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Abstract 

A block of rubber eventually buckles under severe flexure, and sev- 
i— — i eral axial wrinkles appear on the inner curved face of the bent block. 

Experimental measurements reveal that the buckling occurs earlier — 
j>. at lower compressive strains — than expected from theoretical predic- 

ts tions. This paper shows that if rubber is modeled as being bimodular, 

and specifically, as being stiffer in compression than in tension, then 
If} flexure bifurcation happens indeed at lower levels of compressive strain 

than predicted by previous investigations (these included taking into 
O account finite size effects, compressibility effects, and strain-stiffening 

effects.) Here the effect of bimodularity is investigated within the 
• • theory of incremental buckling, and bifurcation equations, numerical 

methods, dispersion curves, and field variations are presented and dis- 
^ cussed. It is also seen that Finite Element Analysis software seems 

to be unable to encompass in a realistic manner the phenomenon of 

bending instability for rubber blocks. 
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1 Introduction 



Rubber is a highly deformable solid. When we are presented with a rub- 
ber band or a rubber string, we almost automatically subject it to a large 
stretch to test its extensibility and then sometimes consider its behavior in 
torsion. However, a thick, stubby, rubber block does not lend itself naturally 
to these two deformations, nor to compression or shear. In fact, one can 
say that large bending (or flexure) is the favorite mode of deformation of 
rubber blocks. Moreover, many engineering and industrial devices rely on 
the excellent bending characteristics of rubber blocks and subject them to 
countless bending/straightening cycles. The importance of bending of blocks 
in applications is the motivation for our study of their stability. 

Gent and Cho (1999) subjected blocks of natural rubber to severe bend- 
ing until they saw creases appearing on the inner bent face. We might expect 
these axial wrinkles to form on this bent face because of our intuitive notion 
of a region of tension and a region of compression in the neighborhood of the 
outer and inner faces of a bent block, respectively. Once circumferential line 
elements contract up to a certain critical stretch A cr , say, the experiments 
of Gent and Cho thus suggest that incremental wavy static deformations 
can be superposed on the primary, non-linear bending deformation, signal- 
ing the onset of instability (in the linearized sense). Gent and Cho (1999) 
assumed that the inner face bending instability should occur at the same crit- 
ical stretch as that of surface instability of an incompressible half-space. For 
the neo-Hookean material, they therefore expected the inner face to buckle 
at A cr = 0.54, and were surprised to see it occur earlier, at A cr = 0.65 ± 0.07. 

Their prediction relied on several assumptions. Several studies have 
tested these assumptions to discover whether modifying any of them would 
significantly increase the value of A cr . We now look at these assumptions 
in turn. First, the incompressibility assumption: does the introduction of 
slight compressibility (as is suitable for natural rubbers) change the critical 
stretch of surface instability? The answer, recently established by Murphy 
and Destrade (2009), is that slight compressibility has little quantitative ef- 
fect on the value of the critical ratio of compression for surface instability; 
if anything, it makes the half-space more stable. Next, the half-space as- 
sumption: does the introduction of a finite size for the block make a big 
difference? The answer is also no, as shown by Haughton (1999) and Coman 
and Destrade (2008), with the critical stretch slightly increased from 0.54 
(half-space in plane strain compression) to 0.56 (block in large bending), 
with very little variation with the block's dimensions. The final assumption 
is that the neo-Hookean strain-energy function is an adequate model of nat- 
ural rubber. In fact, the Mooney-Rivlin material, and also the most general 
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third-order elasticity model of incompressible materials, are all equivalent 
to the neo-Hookean material in the case of plane strain bending (Goriely et 
al., 2008; Destrade et al., submitted) so that quite a wide variety of mate- 
rials reduce to the neo-Hookean form, although it is characterized by only 
one material parameter, the initial shear modulus. Moreover, Gent and Cho 
(1999) are correct in their assertion that it might not be 'necessary to con- 
sider stress-strain relations incorporating finite-extensibility effects' for this 
problem. Destrade et al. (2009a) studied the impact of the strain-stiffening 
effect on bending instability. Their conclusion is that it does promote insta- 
bility for materials which stiffen early (typically at 10-20% extension stretch) 
such as biological soft tissues, but that it does not affect the bending insta- 
bility of materials which stiffen at large stretches (at 200-600%, say) such as 
natural rubbers. 

In this paper, we investigate the influence of a completely different as- 
sumption, one that is usually implicit in the study of non-linearly elastic 
materials and, in particular, implicit in the work of Gent and Cho (1999). 
This assumption is that the mechanical behavior of rubber is the same in 
tension as in compression. 

Intuitively, we expect that there are regions of tension and regions of 
compression in a deformed solid; however it is a far from trivial task to 
provide a rigorous definition of these terms, especially in non-linear elasticity, 
as can be seen in the elegant study by Curnier et al. (1999). Here we focus 
on the flexure of a block, where intuition suggests that the region near the 
outer face of the bent block is a region of "tension", and the region near 
the inner face is a region of "compression". In these regions, line elements 
originally aligned along the length of the bar are extended and contracted, 
respectively. 

Another deformation where the distinction between tension and com- 
pression can be determined on the basis of intuition alone is that of simple 
tension and compression. Subjecting a sample of rubber to a large homoge- 
neous tension is a routine mechanical engineering experiment. Subjecting a 
block of rubber to a large homogeneous compression is, on the other hand, 
a much more difficult task. In particular, it requires generous lubrication to 
avoid non- homogeneous bulging of the sides and this can lead to slippage 
of the block (Brown, 2005) and tilting of the platens. There is a dearth of 
experimental data for both tension and compression tests of a given rubber 
sample. The limited data available in the literature suggest that rubbers can 
behave in a completely different way in one protocol compared to the other. 
For example, Bechir et al. (2006) measure a Poisson ratio of 0.48 in exten- 
sion for NR70 (i.e. close to incompressibility) but a Poisson ratio of 0.26 in 
compression. Such bimodularity has also been observed in a variety of other 
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Figure 1: Bending a block of silicone, until creases appear on the inner face. 
Here the block has dimensions 16.5 cm x 5 cm x 6 cm, giving an aspect ratio 
L/A = 3.36. The end faces of the block were glued onto plexi-glass plates, 
and these are subjected manually to a bending moment. At least three axial 
wrinkles are visible. 
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elastic materials, including rocks (Lyakhovsky et a!., 1997), nacre (Bertoldi 
et al., 2008) and cartilage (Soltz and Ateshian, 2000). A biomechanical ex- 
ample of bimodularity can be found in Mirnajafi et al. (2006) who measured 
the flexural stiffness of pig aortic valve leaflets in the direction of natural 
leaflet motion (Young's modulus E + ) and against that direction (Young's 
modulus E~) and found that 0.43 < E + /E~ < 0.78. Bimodularity, or ten- 
sion/compression asymmetry, has long been recognized as being important 
in the theory of plasticity. Differences in the yield stresses between tension 
and compression have been observed, for example, in copper alloys (Yapici 
et al., 2007; Kuwabara et al., 2009) and in metallic glasses (Schuh and Lund, 
2003). It is our contention that bimodularity can have a significant effect on 
the behavior of elastic materials as well and should be taken into account, 
where possible, in the analysis of their behavior. 

We illustrate our argument with a study of the instability of flexure and 
find that it is promoted when the rubber is stiffer in compression than in 
tension. We take the neo-Hookean model (equivalent to Mooney-Rivlin and 
third-order elasticity models in bending) as the base strain-energy density. 
As mentioned earlier, only one material parameter plays a role in bending 
for these materials, namely the initial shear modulus. Therefore only fx, the 
dimensionless ratio of the shear moduli in tension and in compression, plays 
a role in the determination of the critical stretch of compression of bimodular 
third-order elasticity incompressible solids. We find that the lower and mean 
values of Gent and Cho's experimental range, A cr = 0.58 and A cr = 0.65, are 
reached for a bimodular block with fi = 0.68 and fx = 0.44, respectively. In 
general, we find that bimodularity increases the value of the critical stretch 
ratio of compression, and decreases the number of expected wrinkles. A 
somewhat counterintuitive outcome is that it also increases the allowable 
angle of bending prior to bifurcation. These results (§ |3.3[ ) are the product of 
an in-depth analysis of the large deformation field ($21) and of the incremental 



equations of equilibrium (£3.2). They rely on advanced numerical methods 
and equations developed elsewhere (Destrade et al., 2009a). 

In the concluding section (Q, we also present the output of a Finite 
Element Analysis simulation for the flexure of a uni-modular neo-Hookean 
block, and find that it is not satisfactory when confronted with experimental 
observations (Figure [T]) and theoretical predictions. 

2 Bimodular bending deformation 

Take a block of length 2L, height H, and thickness 2 A: 

-A<X<A, —L <Y<L, 0<Z<H, (2.1) 
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made of an incompressible, isotropic, homogeneous, non-linearly elastic bi- 
modular solid. Bend it under applied terminal moments into the following 
circular, annular sector: 

r a < r < r b , -a<6 < a, < z < H, (2.2) 

where r a , are the inner and outer radii of the bent block's curved faces, 
corresponding to the planes X = —A and X = A in the reference configura- 
tion respectively, and 2a is the bending angle (0 < a < 7r). Here (X, Y, Z) 
and (r, 9, z) are the rectangular and cylindrical polar coordinates of a particle 
before and after deformation, respectively. Flexure is described by Rivlin's 
(1949) solution 



'2LX/a + (r 2 a + r 2 b )/2, 9 = aY/L, z = Z, (2.3) 
from which it follows that 

r\ -r 2 a = 4AL/a. (2.4) 
The circumferential stretch A2 is a useful, non-dimensional quantity, given 

by 

A 2 = ar/L. (2.5) 

If then 

A2 = ar a /L, \\ = an/L, (2.6) 



we may rewrite (2.4) as 
where 



(\*y - (A2) = 2e, (2.7) 



2a{A/L) (2. 



is a non-dimensional measure of the amount of bending. Note that the curve 
separating the region of tension, where 1 < A2 < A2, from the region of 
compression, where \% < A2 < 1 is called the neutral axis, defined by, for 
example, Varga (1966), as 

r = r n = L/a. (2.9) 

We need a second relationship between Ag and \\ in order to fully de- 
termine the deformation field. The classical solution of the flexure problem 
due to Rivlin, where there is of course no bimodularity, assumes that the 
curved surfaces are stress-free and this yields A2A2 = 1, which is independent 
of the strain-energy function. In a bimodular block, the stress-free boundary 
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condition yields a second relationship that does depend on the form of the 
strain-energy function assumed. Here we assume that 

H^ = ^ (A 2 2 + A 2 - 2) , (2.10) 

where ^ are the shear moduli in tension and compression respectively. This 
is the plane strain form of the neo-Hookean strain-energy function and we 
recall from the Introduction that a wide class of materials reduces to this form 
in plane strain. For this material, Destrade et al. (2009b) show that when 
the curved faces of the bent block are free of normal traction and the stress 
is continuous throughout the bimodular block, then the following equation 
is obtained: 

2ex 3 + [l - 4e - /2 (1 - 2e) 2 ] x 2 - 2 [1 - e - jl (1 - 2e)] x + 1 - /2 = 0, (2.11) 

where x = (A 2 )~~ 2 and fi = fi + /fi~ is the ratio of the shear moduli in the 
regions of tension and compression, respectively. This is a bicubic in A 2 and 
a quadratic in e. 



Equations (2.7) and (2.11) determine completely the large bending defor- 



mation for a bimodular neo-Hookean block, with a given fi, once the amount 



of bending e = 2a(A/L) is prescribed, since A 2 can be determined from (2.11 ) 
and A 2 then obtained from (2.7). To illustrate the magnitude of the stretches 



obtained, Figure [2] shows the variations of A 2 with e. 

Once e, Af , and A 2 are known, the inner and outer radii, as well as the 
neutral axis, are completely determined. Normalized with respect to the 
block thickness, they are 

r^_A| _ 1 J±_^2 / 919 x 

2A ~ e ' 2A ~ e' 2A ~ e ' [ ' 

Finally, the radial stress field is either [i + a + or [i~a~ , depending on the 
region — tension or compression. Here, the non-dimensional quantities <j ± are 

a+ = \ [A 2 + A 2 2 - (A2) 2 - (A fe 2 )- 2 ] , a = \ [A 2 + A 2 2 - (A 2 ) 2 - (A^)" 2 ] 

(2.13) 

This field leaves the bent faces free of normal traction. 



3 Bimodular bending instability 



3.1 The general method 



Suppose that instead of specifying e, as in the last section, A 
ential stretch on the inner face, is prescribed. Then (2.11) can be solved to 



2, the circumfer- 
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Figure 2: Variations of Ag, the circumferential stretch ratio on the inner bent 
face, with e, the product of the block aspect ratio and the bending angle, 
for bimodular solids which are stiffer in compression than in tension. The 
ratio of the tension shear modulus to the compression shear modulus takes 
the values: fi = 1.0 (dotted line; homogeneous block), 0.8, 0.68, 0.6, 0.5, 
and 0.44 (full lines). The circles give an indication of the critical stretch of 
bending instability, see Figure [3] 
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determine e, with A 2 determined from (2.7). We can then investigate whether, 
at the prescribed stretch A 2 , a static incremental field can be superimposed 
onto the large bending solution. If so, then A 2 is what is termed a critical 
stretch of bending instability, A cr . The incremental equations of equilibrium 
used to determine these critical stretches are described next. 



3.2 Bifurcation 

For incompressible solids, the incremental equations of equilibrium are (Og- 
den, 1984) 

divs = 0, divu = 0, (3.1) 

where s is the incremental nominal stress and u is the incremental mechanical 
displacement. We seek solutions of these equations in the form 

{u, s} = ${ [U(\ 2 ), V(\ 2 ), 0, S rr (\ 2 ), S re (\ 2 ), 0] e ine } , (3.2) 

where U, V, S rr , S r g, are complex functions of A 2 only and 

n = pn/a, p an integer, (3.3) 

gives the number of wrinkles on the bent faces. 

Destrade et al. (2009a) formulated the following Stroh form of the incre- 
mental equations: 

^-T7(A 2 ) = ^G±(A 2 )r7(A 2 ), (3.4) 

where 

f 1 =[U,V,irS rr: irS re \ t (3.5) 

is the displacement-traction vector (Shuvalov, 2003). The components of the 
matrices are given by 



where 



G 



± 



—n 



-n(l-a ± \ 2 2 ) -i(l- ( x ± A|) 
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(3.7) 
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The bifurcation problem is a two-point boundary value problem. Thus 



(3.4) must be integrated numerically, subject to the following boundary con- 



ditions of zero incremental traction on the bent faces: 

= [U(X a 2 ), Of, ry(A^) = [U(X b 2 ), 0]\ (3.8) 



Following the approach of Haughton (1999), the value of n, defined in (3.3), 
is determined by assuming that there are no incremental normal tractions on 
the plane end-faces 9 = ±<x We must also enforce the continuity of displace- 
ment and traction (i.e. of r/) across the neutral axis (at A2 = 1). Several 
strategies exist to tackle this problem. To find A cr only, the compound matrix 
method is most efficient; to find A cr and the mechanical fields throughout the 
block, the impedance matrix method works very well. In the former case, 
we have to integrate numerically a linear system satisfied by the compound 
matrix; in the latter case, we have to integrate numerically a Riccati equa- 
tion satisfied by the conditional surface impedance matrix (see Destrade et 
al. (2009a) for details). 

3.3 Results 

For a given bimodular material, fi = fi + /fi~ is prescribed. Non-dimensionali- 



sing the incremental equations (3.4) yields two independent, non-dimensional 
parameters: A cr and L/(pA). Plotting then A cr versus L/(pA) yields master 
curves for the "physical" dispersion curves, which give the critical stretch 
ratio A cr versus the block aspect ratio L/A, for p — 1, 2, 3, . . .. The maximum 
of the master curve is a good indicator of the physical critical stretch of a 
block with a given aspect ratio. Further details can be found in Haughton 
(1999), Coman and Destrade (2008), and Destrade et al. (2009a). 

Consequently, we plot these master curves for various values of fi within 
the range 0.44 < p, < 1; the upper value restricts attention to bimodular 
solids which are stiffer in compression than in tension, the choice of the 
arbitrary lower value will be explained shortly. Qualitatively, we find that 
the bending instability occurs earlier for bimodular solids than for a uni- 
modular block (/a = 1), in line with results on coated blocks (Dryburgh 
and Ogden, 1999) and bilayered blocks (Roccabianca and Bigone, private 
communication). Quantitatively, we find that the lower value of the Gent 
and Cho (1999) range, i.e. A cr = 0.65 — 0.07 = 0.58, corresponds to fi = 0.68 
whilst their mean value, A cr = 0.65 corresponds to ft = 0.44, see Figure |3j 
On that figure, we note that the maximum of each master curve increases as 
the bimodularity becomes more pronounced. 

There is a seemingly counter-intuitive result worth noting: although insta- 
bility is clearly promoted by bimodularity, in the sense that A cr is increased, a 
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0.5 1.0 1.5 2.0 2.5 3.0 

U(pA) 

Figure 3: How bimodularity promotes bending instability. As the differ- 
ence between compressive and tensile behaviors becomes more pronounced 
(/} = 1.0, 0.8, . . . , 0.44), the critical stretch of bending increases. (We checked 
separately that when the solid is stiffer in tension than in compression, the 
block becomes more stable.) The circles indicate the maximum of each curve, 
and give a good indication of the "physical critical stretch" of bifurcation, 
see Figure |4| 
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bimodular block can nonetheless be bent further than its uni-modular coun- 
terpart before instability occurs, given that instability for a bimodular block 
corresponds to a larger value of the amount of bending e = 2a(A/L). This 
is clear from the observation of the circles on Figure |2j which indicate the 
onset of bending instability. This occurs not alone at a higher value of A2 
than the homogeneous (uni-modular) block (at fx = 1), but also at a higher 
value of e. We illustrate this observation further with a specific example. 

Figure [4] displays the dispersion curves p = 1, 2, . . . , 6 when fx = 0.68. 
When the aspect ratio L/A = 3.36 (this is the aspect ratio of the silicone 
block used in Figure [T]), it is clear from the figure that its inner face buckles 
in bending when p = 3 and A2 = 0.5797. 

Figure [5] contrasts this incremental buckling field throughout the bent 
block (right-hand side of the figure) with the corresponding uni-modular 
field. From (\2A1 ), the corresponding amount of bending for the bimodular 
block is e = 1.660, giving a bending angle of 2a = e{L/A) = 319°. In 
comparison, a uni-modular block {ft = 1) with the same dimensions buckles 
in bending with p = 8 at A cr = 0.5612 (left side of the figure). Hence, 
by moving the maximum of the master curve upwards and to the right on 
Figure [3j bimodularity promotes earlier buckling (A cr increases), with fewer 
wrinkles (p decreases). At the same time, e increases too (see Figure ^), 
which means that the critical bending angle is larger. Note that Figure 5 
displays the correct (scaled) size of the bent blocks relative to the undeformed 
blocks. Similarly, the amplitude of the incremental wrinkles on the outer 
face is computed relative to their amplitude on the inner face. However, the 
amplitude of the wrinkles relative to the large bending displacement is not 
known, because the incremental analysis is linear. 



4 Concluding remarks 

Gent and Cho (1999) showed experimentally that a block of rubber buckles 
in bending when the circumferential stretch is 0.65 ± 0.07, earlier than their 
prediction based on surface buckling, where the critical stretch is 0.54. Pre- 
vious modeling attempts at finding a satisfactory explanation for this early 
buckling have failed. They included taking into account slight compressibil- 
ity, finite dimensions, and the strain-stiffening effect. Here we investigated 
how bimodularity might affect the bending instability of rubber blocks, and 
found that if the rubber is stiffer in compression than in tension, then the 
bending buckling does occur earlier. Our detailed analysis also revealed that 
the bimodular block can be bent by a larger angle than the corresponding 
uni-modular block with the same dimensions. It was important to evaluate 
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Figure 4: Dispersion curves for a given bimodular block, with fi = 0.68. 
Depending on the aspect ratio L/A of the block, a certain mode of buckling 
is triggered at A cr . Except for unrealistically short and stubby blocks, A cr is 
confined within the range 0.578 — 0.580, that is, close to the value indicated 
by the circle on Figure [3] For example, a block with aspect ratio L/A = 3.36 
(like the silicone block of Figure [I]) experiences bending instability when 
the circumferential stretch reaches 0.5797 approx. and three axial wrinkles 
should appear (p — 3), as can be easily seen from consideration of the vertical 
hash line. 
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Figure 5: Mechanical displacement field in bent blocks with aspect ratio 
L/A = 3.36 (like the block of Figure [I]). On the left: uni-modular block 
(jl = 1); on the right: a bimodular block with fi = 0.68. First stage: unde- 
formed configuration. Second stage: maximum bending, prior to bifurcation; 
here the angles of bending are 275° and 319°, for the uni-modular and the 
bimodular blocks, respectively, corresponding to critical stretches of 0.56 and 
0.58, respectively. Third stage: incremental buckling; here 8 wrinkles and 
3 wrinkles appear on the inner face of the uni-modular block and of the bi- 
modular block, respectively. The thick line is the neutral axis, where line 
elements are neither extended nor contracted. 
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the effect of bimodularity for the flexure stability problem, because numeri- 
cal simulations using the Finite Element Method do not provide physical and 
realistic predictions. We conclude the paper by summarizing our numerical 
simulation experiments. 

We used the commercial Finite Element Analysis software Abaqus [Das- 
sault Systemes Simulia Corp., Providence, Rhode Island, USA] to simulate 
numerically the flexure of a uni-modular block with the same dimensions 
(16.5 cm x 5 cm x 6 cm) as the silicone block of Figure 1. We measured its 
mass density as being p = 641.176 kg.m -3 . A hyperelastic material model 
was chosen using the neo-Hookean form of the strain energy potential defined 
as follows: initial shear modulus fi = 0.6 MPa and initial Poisson ratio 0.4999 
(to mimic incompressibility) . A reference point at the center of each of the 
respective end faces of the block governed the movement of each of these 
faces by means of a coupling constraint. Each of the reference points was 
displaced 3.91 cm towards the other, and they were rotated by 60°. These 
values were computed as being consistent with Rivlin's plane strain solution 
(2.3). The displacement and rotation were ramped linearly over a time step 
of 0.5 s. The mesh consisted of first-order fully integrated 8 noded hexahedral 
3D continuum elements, enhanced by incompatible modes to improve their 
bending behavior. The analysis was computed using explicit procedures to 
allow for the large deformation. We looked at different types of element for- 
mulations and mesh refinements, but always found the same following trends 
in the final output. As shown on Figure |6j Abaqus predicts that under severe 
bending, the block should experience the formation of a single fold, centered 
axially on the inner bent face. This fold occurs at a smaller bending angle, 
and at larger strains, than observed experimentally and predicted theoreti- 
cally. Other inconsistencies include the absence of other folds, the opposite 
of the anti-clastic effect (in the sense that the outer face widens while the 
inner face shrinks), and the penetration of the fold inside the block (as seen 
on the last picture). 
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Figure 6: Finite element simulation of the large flexure of a block with the 
same dimensions as the block of Figure [TJ The top pictures shows the results 
of the simulation just prior to folding of the inner surface (first picture) and 
after (second picture). The bottom pictures show zooms on the fold in the 
buckled geometry. 
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